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ABSTRACT 

Context. The atmosphere of hot Jupiters can be probed by primary transit and secondary eclipse spectroscopy. Due to the intense UV 
irradiation, mixing and circulation, their chemical composition is maintained out of equilibrium and must be modeled with kinetic 
models. 

Aims. Our purpose is to release a chemical network, and the associated rate coefficients, developed for the temperature and pressure 
range relevant to hot Jupiters atmospheres. Using this network, we study the vertical atmospheric composition of the two hot Jupiters 
(HD 209458b and HD 189733b) with a model that includes photolyses and vertical mixing and we produce synthetic spectra. 
Methods. The chemical scheme is derived from applied combustion models that have been methodically validated over a range of 
temperatures and pressures typical of the atmospheric layers influencing the observations of hot Jupiters. We compare the predictions 
obtained from this scheme with equilibrium calculations, with different schemes available in the literature that contain N-bearing 
species and with previously published photochemical models. 

Results. Compared to other chemical schemes that were not subjected to the same systematic validation, we find significant differences 
whenever non-equilibrium processes take place (photodissociations or vertical mixing). The deviations from the equilibrium, and thus 
the sensitivity to the network, are more important for HD 189733b, as we assume a cooler atmosphere than for HD 209458b. We found 
that the abundances of NH 3 and HCN can vary by two orders of magnitude depending on the network, demonstrating the importance 
of comprehensive experimental validation. A spectral feature of NH 3 at 10.5 yum is sensitive to these abundance variations and thus to 
the chemical scheme. 

Conclusions. Due to the influence of the kinetics, we recommend the use of a validated scheme to model the chemistry of exoplanet 
atmospheres. The network we release is robust for temperatures within 300-2500 K and pressures from 10 mbar up to a few hundreds 
of bars, for species made of C, H, O, N. It is validated for species up to 2 carbon atoms and for the main nitrogen species (NH 3 , HCN, 
No, NOx). Although the influence of the kinetic scheme on the hot Jupiters spectra remains within the current observational error 
bars (with the exception of NH 3 ), it will become more important for atmospheres that are cooler or subjected to higher UV fluxes, 
departing more from equilibrium. 

Key words. Astrochemistry - Planets and satellites: atmospheres - Planets and satellites: individual: (HD 209458b - HD 189733b) 



1. Introduction 

So far, more than 700 exoplanets have been confirmed and thou- 
sands of transiting candidates have b een identified by the space 
telescope Kepler dBatalha et al.l2Q12f) . Among them, hot Jupiters 
are a class of gas giants with orbital periods of a few days or 
less. They a re found around ~ 0.5% of KGF stars dHoward et alj 
l2rMl2012h . About 10% of them transit their host star and their 



ied by transit spectroscopy (e.g. Charbonneau et al. 2000, 2008 


Richardson et al. 120071 Tinetti et all 2007 


Sing et al. 


2008 


Swain et al.ll2008allbl 2009allbl Huitson et al. 


2012). 



Although current observations are still limited and sub- 
jected to diver gent interpretations, future i nstruments such a s 
E-ELT, JW ST fcardner et al.ll2006h . EChO (iTinettiet al.ll201lh . 
FINESSE dSwainl 12010) should provide better constraints on 
both the chemical composition and the temperature profiles of 
the nearby hot Jupiters like HD 189733b and HD 209458b. They 



will also be able to study more distant targets and deliver statis- 
tically significant trends about the nature of their atmospheres. 
Chemical modeling will be an important component of these 
studies. It will point to key observations able to distinguish be- 
tween various hypotheses and will be used to analyze the obser- 
vations and constrain, for instance, the atmospheric elemental 
abundances. 

The first models of hot Jupiter atmospheres assumed chemi- 
cal equilibrium (e.g.lBurrows & Shartir i999: Seager & Sasselov 
2000; Sharp & Burrows! l2007t iBarmanl 120071: iBurrows et al. 
20071 120081: iFortnev et all l2008al) . However, strongly irradi- 
ated atmospheres are unlikely to be at chemical equilib- 
rium. Their intense UV irradiation (typically 10,000 times the 
flux received on the top of the atmosphere of Jupiter) and 
strong dynamics result in photolyses and diffusion/advection 
timescales that are comparable or shorter than the chemical 
ones. Deviations from the thermodynamic eq uilibrium have been 
discussed with timescale arguments (e.g. Fodders & Feglev 
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2002; For tnev et al.l l2006t l2008ht IVisscher et all 120061 120101: 
Madh usudhan & Seagerll2010l) . or modeled with a few reactions 
describing the CO-CIT, co nversion coupled with the dynamics 
dCooper & Sho wman 2006). A more detailed modeling requires 
the use of a photochemical kinetic network. A kinetic network 
is, in practice, a list of reactions and associated rate coefficients 
able to describe quantitatively (within a certain accuracy) the 
kinetics of a pool of species, usually the most abundant ones. 
Constructing such network of reactions implies to solve two ma- 
jor questions. One has to do with the completeness of the net- 
work: What are the species and the reactions connecting them 
that must be included? The other issue is the availability of the 
kinetic data: The literature and databases may not provide the 
rate coefficients for some of the needed reactions or may pro- 
vide conflicting values with no recommendation. These two is- 
sues are tightly connected and both depend on the considered 
range of temperatures and pressures. Eventually, and whatever 
the methodology adopted to select the reactions and their rates, 
it is the ability to predict experimentaly-controled abundances 
that can validate or not the network. 

To investigate the consequences of the strong UV incident flux 
on ne utral species, photochem ical models have then been devel- 
oped dLianget al.l 120031 [2004). Based on kinetics model ded- 
icated to Jupiter's low-temperature atmosphere, these models 
however neglected endothermic reactions, which a re in fact 
rather efficient in such hot atmospheres. iLine et all d2010l) in- 
troduced some endothermic reactions to a similar Jovian pho- 
tochemical scheme but most of the pre-existing reactions were 
not reversed. They were therefore not able to reproduce the ther- 
modynamic eq uilibrium, wh i ch occu rs in the deep atmospheres 
of hot Jupiters. IZahnle et al.1 d2009allbl) developed a photochem- 
ical model considering the reversal of their whole set of two- 
bodies exothermic reactions. They selected their rate constants 
in the NIST database^ based on the following criteria: relevance 
of temperature conditions, date of review, date of the experi- 
ment and da t e of th e theoretical study (in order of preference). 
Mos es etall d201 ll) developed a model that considered the re- 
versal of all the reactions, including three-body reactions, ensur- 
ing to reproduce the thermodynamical equilibrium from the top 
to the deepest layers of the atmosphere. Their c hemical scheme 
is derived from the Jupiter and Saturn m odels ([Gladstone et alj 
1 1 9961: iMosesll 1 9961: iMoses et al Jl 1 995bllal l2000albT) " with further 
updates on the basis of combustion-chemistry literature. None of 
the aforementioned works discuss the validation of the chemical 
scheme against experiments. In addition, the fact that computed 
abundances evolve towards the composition predicted by equi- 
librium calculations (at given pressure P and temperature T with 
no external irradiation nor mixing) is by no means a validation 
of the kinetic network. Indeed, any network containing at least 
as many independent reversible reactions as modeled species, 
in which the rates for the backward processes are derived from 
equilibrium constants and forward rates, will evolve toward 
the equilibrium predicted with the same equilibrium constants, 
whatever the quantitative values of the forward rates, as illus- 
trated in Fig.Q] Fortunately, and due to the physical conditions 
and elemental composition of hot Jupiters (and hot Neptunes) 
atmospheres, we benefit from decades of intensive work done in 
the field of combustion, that includes a vast amount of experi- 
ments, the development of comprehensive mechanism^ and the 




Fig. 1: Abundances of NH 3 as a function of time computed with 
two kinetic schemes fully reversed according to equilibrium con- 
stants but differing by their nitrogen chemistry (nominal and 
GRIMECH as defined in Sects. 12X21 and l33l While they both 
converge towards the equilibrium (dotted line) they exhibit very 
different evolution. Initial condition is a mixture of H 2 , CH 4 , 2 , 
N 2 and He with solar elemental abundances. 



systematic comparison between the two. Therefore, we propose 
in the present work a new mechanism dedicated to the chemical 
modeling of hot atmospheres that is not adapted from previous 
Solar System photochemical models but derives from industrial 
applications (mainly combustion in car engines). Details about 
this chemical network and its range of validity are presented in 
Sect.O 

We use this chemical network in a ID model that includes 
photolyses and vertical transport, which has been previously 
used to study the atmospheric photoc hemistry of various ob - 
jects of the Solar System: Ne ptune dDobriievic et al.l 1201 0l) . 
Titan dHebrard et al.l2006ll2007l) . Saturn dDobriievic et alj 2003: 
ICavahe et all l2009l) . and Jupiter (ICavahe et all 120081) ') as well 
as extrasolar terrestrial planets dSelsis et al.l 120021) . We model 
the photochemistry of two hot Jupiters: HD 209458b and 
HD 189733b (Sect. I231 l, We study the departure from ther- 
modynamic equilibrium and compare our results with those of 
Mos es et alj(l2011h (Sect.[3]l. We also investigate how including 
different reaction networks specific to nitrogen-bearing species 
influences the model results (Sect. 13. 3t and the planetary syn- 
thetic spectra (Sect. [3T3~Tl >. 

2. The model 

2. 1 . Kinetic network: from car engine to hot Jupiters 

Significant progress has been done during the past decade in 
the development of validated combustion mechanisms. In the 
context of limiting the environmental impact of transportation, 
there is indeed a need in the development of detailed chemi- 
cal kinetic models more predictive and more accurate for the 
combustion of fuels. One part of the studies undertaken in the 
LRGP (Laboratoire Reactions et Genie des Procedes, Nancy, 
France) concerns engine-fuel adaptation in order to improve the 
efficiency of engines and to limit the emission of pollutants. 
Gasoline and Diesel fuels contain a large number of molecules 
belonging to several major hydrocarbon families. Biofuels con- 
tain also oxygenated species, such as alcohols and methyl es- 



1 http://www.kinetics.nist.gov/kinetics 

2 In the field of combustion, a mechanism or reaction base is a net- 
work of reactions able to describe the kinetic evolution of a given pool 
of species. The mechanism includes the list of reactions and the associ- 



ated rate coefficients, in a modified Arrhenius form, as well as the ther- 
modynamic data for all the species involved in these reactions, which 
are required to calculate the equilibrium constants of the reactions and 
the rates of the reverse reactions 
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ters. Oxidation and combustion of these complex blends occur 
by radical chain reactions involving hundreds of molecular and 
radical species and several thousands of elementary reactions in 
the case of pure reference fuels such as n-heptane, iso-octane 
or cetane. The primary focus of the currently developed chemi- 
cal models is to simulate the main combustion parameters (auto- 
ignition delay times, laminar flame speed, heat release), which 
are needed for the design of engines or burners, to estimate the 
fuel consumption, and to model the formation of some of the 
main regulated pollutants (carbon monoxide, nitrogen oxides, 
unburned hydrocarbons and particulate matter). Most of these 
kinetic models were developed for industrial applications and 
have been validated in a range of temperatures, from 300 to ap- 
proximately 2500 K, and for pressures from 0.01 bar to some 
hundreds of bar. What is worth noticing is the similarity of these 
temperature and pressure ranges with the conditions prevailing 
in hot Jupiters atmospheres, in the very layers where they in- 
fluence the observed molecular features. In addition, combus- 
tion mechanisms mainly deal with molecules made of C, H, O 
and N, which are also the main constituents of the molecules 
and radicals in these atmospheres. For this reason, we have de- 
cided to implement such a mechanism, which has already been 
appli ed successfully to many cases and systematically validated 
dBounaceu r et al. 2007), to study the atmosphere of hot Jupiters. 

In this study we have used a C/H/O/N mechanism, whose 
core is a C0-C2 mechanism that includes all the reactions re- 
quired to model the kinetic evolution of radicals and molecules 
containing less than three carbon atoms. This mechanism also 
contains some species with more carbon atoms that are needed 
to model the abundance of C0-C2 species. This mechanism does 
not include nitrogen species, except N 2 as a third body. As ni- 
trogen species, such as N 2 , NH 3 , HCN, CN, are expected to 
be important constituents of hot Jupiter atmospheres, we com- 
pleted this C0-C2 base with a validated sub-mechanism specif- 
ically constructed to model nitrogen species and all the cross- 
term reactions involved (for instance, reactions between alka- 
nes and NOx). These mechanisms do not use rate coefficients 
that have been adjusted by optimization procedures in order to 
fit experiments. Their values are those recommended for the in- 
divid ual processes by the main kinetics databases for combus- 
tion (iTsang & Hampson|[T986t iManion et al.ll2008t [Smith et al. 
ll999tlBaulch et al.1120051) . The list of the reactions and their rate 
coefficients are available in the online database KID A: Kinetic 
Database for Astrochemistr>0 (Wakel am et al.ll2012h . The final 
mechanism includes 957 reversible and 6 irreversible reactions 
(see Sect. 12. 1.3l , involving 105 neutral species (molecule or rad- 
ical). Helium is also included in this mechanism and plays the 
role of third body in some reactions. 

2.1.1. C0-C2 reaction base 

The C0-C2 reaction base we use was develo ped for in- 
dustrial applications and was first presented by Barbe et al.l 
( 1995b and has been co ntinuously updated dFournet et al.| [1999; 
Bounaceur et al . 2010). This mechanism is designed to repro- 
duce the kinetics of species with less than three carbons. It in- 
cludes all the unimolecular or bimolecular reactions involving 
radicals or molecules containing no more than two carbon atoms. 
This mec hanism has been built by using a reaction grid, as pro- 
posed by Tsang & Hampson ( 1986). Every unimolecular and bi- 
molecular elementary forward reactions involving the consid- 
ered reacting species have been systematically written. Reacting 
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species include 46 compounds (19 molecules and 27 radicals) 
which have been ranked according to the molecular formula 
OtCjH, (with x varying from to 3, y from to 2 and z from 
to 6): CO, H 2 , H 2 0, 2 , H 2 2 , CH 4 , H 2 CO, CH 3 OH, C0 2 , 
CH3OOH, C 2 H 2 , C 2 H 4 , C 2 H 6 , CH 2 CO, CH 3 CHO, C 2 H 5 OH, 
C 2 H 5 OOH, CH 3 COOOH, cC 2 H 4 (Ethylene Oxide), C, CH, 
'CH 2 (singlet), 3 CH 2 (triplet), 0( 3 P), H, OH, OOH, CH 3 , HCO, 
CH 2 OH, CH 3 0, CH 3 00, C 2 H, C 2 H 3 , C 2 H 5 , CHCO, CH 2 CHO, 
CH 3 CO, C 2 H s O, C 2 H 4 OOH, C 2 H 5 00, CH 3 COOO, CH 3 OCO, 
C0 2 H, l-C 2 H 4 OH and 2-C 2 H 4 OH (ethyl radical isomers, 1- 
hydroxy and 2-hydroxy). This ranking permits to separate easily 
the part of the mechanism related to pyrolysis reactions from 
the one related to oxidation or combustion. The mechanism 
also includes 14 species containing 3 or 4 carbon atoms: C 3 H 8 , 
C 4 H 8 , C 4 H 10 , C 2 H 5 CHO, C 3 H 7 OH, C 3 H 7 0, C 4 H 9 0, C 2 H 6 CO, 
C 3 H 8 CO, C 2 H 3 CHO, n-C 3 H 7 , i-C 3 H 7 (isopropyl and n-propyl 
radical isomers), 1-C 4 H 9 and 2-C 4 H 9 (1 -butyl and 2-butyl radi- 
cal isomers). 

This C0-C2 mechanism has been widely validated in 
the 300-2500 K, 0.01-100 bar range for several types 
of reactors such as shock-tubes, perfectly stirred reac- 
tors, plug-fl ow reactors, rapid compre s sion machines, laminar 
flames (e.g. iBattin-Leclerc et al.l 12006]: iBounaceur et all 1 2007 ; 
lAnderlohr et al.ll2010tlBounaceur et al.l2010HWang et alj2010l) . 
Obviously, it is not possible to describe in details all these valida- 
tions, but we can mentio n, for instance, the very recent work of 
iDirre nberger e t al.l (1201 ll) who has studied experimentally and 
modeled with success the laminar burning velocity of several 
mixtures including air, hydrogen and components of natural gas. 
Laminar burning velocities are important parameters in many ar- 
eas of combustion science such as the design of burners and the 
prediction of explosions. They also play an essential role in de- 
termining several important aspects of the combustion process 
in spark ignition engines. These experiments have been done 
in specific mixtures, containing nitrogen in the sole form of N 7 
and in which nitrogen species produced from N 2 (NOx in the 
typical mixtures used in combustion) do not significantly affect 
the outcome, in order to validate the C0-C2 mechanism itself. 
Therefore, a model including only the C0-C2 base would not be 
accurate to predict the abundance of C0-C2 species in this range 
of P and T when applied to mixtures containing or producing 
(by reaction with N 2 ) significant levels of nitrogen species other 
than N 2 . 



2.1.2. Nitrogen reaction base 

In our nominal model, the sub-ne twork f o r the nitro- 
gen bearing specie s is d erived from iKonnovl (|2000, 2009) 
and [Coppens e t al.l (I2007I) . It is based on a comprehensive 
analysis of the combustion chemistry of nitrogen oxides 
(IKonnov & De Ruvckl 1 1999a), ammon ia ( Konn ov & De Ruvckl 
2000b), hydrazine (Ko nnov & De Ruvckll200 1b). and modeling 
of nitrogen oxides formation in diff e rent c ombustion systems 
(IKonnov & De Ruvckl fl999b. 2000ai l2001ah . The mechanism 
was tested at the California Institute of Technology, USA, and 
found suitable for steady one-dime nsional detonation and con - 
stant volume explosion simulations (Schultz & Shepherd 1999). 
It was also preferred by the researchers at the University of 
Bielefeld, Ger many to analyze flame st ructure and NO reburning 
in C3 flames (lAtakan & Hartl ieb 2000). In addition ; we consider 
a few additional pathways for HCN oxidation from Dagaut et al. 
(2008). 
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Validations of our nominal sub-network for nitrogen bearing 
species have been made on the basis of experimental data ob- 
tained, for instance, by oxidation of HCN in a silica jet-stirred 
reactor (JSR) at atm ospheric pressure and from 1000 to 1400 K 
(Dag aut et al.ll2008|). or studying lamin ar flame speeds in NH 3 
- N 2 mixtures ( Brown & Smith! [l 9941) . The nitrogen mecha- 
nism includes 42 species (molecule or radical): N0 3 , HON0 2 , 
CH 3 ONO, CH 3 N0 2 , HN0 2 , CH 3 NO, N0 2 , HONO, HCNN, 
HCNO, N 2 0, NCO, HNO, HOCN, NNH, H 2 CN, N( 4 S), CN, 
HNCO, NO, NH, NH 2 , HCN, NH 3 , N 2 , N 2 4 , N 2 3 , N 2 H 2 , 
N 2 H 3 , N 2 H 4 , HNNO, HNOH, HN0 3 , NH 2 OH, H 2 NO, CNN, 
H 2 CNO, C 2 N 2 , HCNH, HNC, HON and NCN. 

For comparison, we have also used other nitrogen sub- 
mechanisms, which are presented in Sec. 13.31 with the corre- 
sponding results. 

Because the mechanism we use was created from individ- 
ual processes and validated without any optimization of their 
reaction coefficients, its application outside the condition range 
of validation is not problematic. This is an issue, for instance, 
with the well-kno wn combustion mechanism GRI-Mech 3.C0 
dSmith et alj 19991) . proposed by Gas Research Institute, which is 
an optimized mechanism designed to model natural gas combus- 
tion. Optimization makes the model extremely accurate within 
the optimization domain but its application beyond is risky 
(iBattin-Leclerc et alJl201 ll) . 



and of hydrogen peroxide into two hydroxyl radicals: 



H 2 2 



OH + OH 



(3) 



At high temperature, chemical kinetics is very sensitive 
to these three reactions, which have been widely studied ex - 
perimentally dBaulch et all [19941: iGoldenl l200l iTroel [20TTh . 
Therefore, we use the kinetic data measured experimentally for 
the forward and the reverse directions instead of calculating the 
reverse rate constants using thermodynamic parameters. 



2.1.4. Excitation of oxygen and nitrogen atoms 

Photodissociations produce excited states of oxygen (O('D)) 
and nitrogen (N( 2 D)) that are not treated in the original com- 
bustion mechanisms. Therefore, we added to the C/H/O/N 
mechanism, 19 reversible reactions which describe the ki- 
netics of O('D) and N( 2 D), including radiative and col- 
lisional desexcitation. These reacti o ns rat es are ta k en (o r 
have been esti m ated) fr om lOkabd dl978l): iHerronl (fl 999): 



lUmemoto et all d!998l): iBalucani et all d2000al): ISatoetal.1 



( 1999); Baluca ni et alj fcoOObl) and lSander et al.1 d201 ll) 



2.1.3. Reversible reactions: kinetics vs thermodynamics 

For most reversible reactions, rate coefficients are only available 
for the exothermic (forward) direction. The rate constant for the 
endothermic (reverse) direction, k r (T) is then calculated as the 
ratio between the forward rate constant kf(T) and the equilib- 
rium constant K eq {T), calculated with thermodynamic al data, as 
explained in AppendixlAl However, rate coefficients have some- 
times been measured for both directions. In such cases, the ratio 
kf(T)/k r (T) departs from K eq (T) as different uncertainties affect 
the rate coefficients and the thermodynamic data. The computa- 
tion of k r (T) using K eq {T) ensures the consistency between ki- 
netics and thermodynamics by making the kinetic model evolve 
strictly toward the thermodynamic equilibrium that we calculate. 
Nevertheless, this choice may not always be the best. It results 
in the propagation into k r (T) of both the errors affecting kf(T) 
and K eq (T). Indeed, thermodynamic parameters used to calcu- 
late K eq (T) are not free of error, and are regularly updated just as 
kinetic data. In the field of combustion, for small species as CH 4 , 
CH 3 and OH, it is common to use experimentally measured ki- 
netic rates, rather than thermodynamical reversal, when they are 
available in the relevant temperature range. There is no obvious 
rule in this matter, but validation of the mechanism with out-of- 
equilibrium experiments seems the only practical way to chose 
between different rates. This is the criterion that we use and our 
nominal network uses thermodynamical reversal for most of the 
reactions but not for three important ones. These reactions affect 
the unimolecular initiations (or thermal dissociation reactions) 
of methane into methyl and hydrogen radicals: 



CH, 



M 



CH 3 + H 



of ethane into two methyl radicals: 



C 2 H 6 



M 



CH 3 + CH 3 



(1) 



(2) 



2.2. Test on the chemical scheme with a 0D model 

In addition to our ID model, we have also developed a simple 
0D model that computes the chemical evolution of a mixture 
at constant temperature and pressure. It does not include mix- 
ing with another mixture nor photolyses. We have used this 0D 
model to compare the composition found at steady state with 
the abundances at thermodynamic equilibrium (calculated with 
the code TECA, described in Appendix|B]i for several couples of 
pressure-temperature. Fig.|2]illustrates this comparison with four 
species. First, we used a version of our nominal scheme in which 
all the reactions are reversed in agreement with their equilib- 
rium constant (solid lines). The computed abundances converge 
exactly towards the equilibrium values (dotted lines) with negli- 
gible numerical differences. Then, we used our nominal model 
in which some reactions are not reversed according to their equi- 
librium constant but using rate coefficients measured experimen- 
tally . In this case, the abundances reached at steady state (dashed 
lines) departs from the predicted equilibrium. This departure re- 
mains very small: below 1% for most species, and always below 
5%. However, one can see that the kinetic evolution can be sig- 
nificantly different, both in terms of abundances and timescales. 

2.3. The 1D model 

To model the chemical composition of the atmosphere of the hot 
Jupiters HD 209458b and H P 189733b, we use our ID time- 
dependent model described in lDobriievic et al.l (12010). As an in- 
put of the model, we give a pressure-temperature profile for the 
atmosphere of the planet studied. This profile is then divided in 
discrete uniform layers with a thickness Az = ^P, where H(z) is 
the pressure scale height. The grid contains ~ 300 layers. Then 
the ID kinetic model resolves the continuity equation (Eq. |4ji as 
a function of time, for each species and atmospheric layer, until 
a steady state is reached. 



http://www.me.berkeley.edu/gri_mech/ 



dn i _ P 
dt ' 



(4) 
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Fig. 2: Comparison between the thermodynamic equilibrium (dotted line) and the evolution of some molecular abundances in the 
OD model, with two different schemes: the thermochemically reversed model (full line) and the nominal model (dashed line) as a 
function of integration time at different temperature-pressure points. Initial condition is a mixture of H 2 , CH 4 , 2 , N 2 and He with 
solar elemental abundances. 



where n, the number density of the species i (cm ), P, its pro- 
duction rate (cm 4 .s _1 ), L,- its loss rate (s _1 ), and O, its vertical 
flux (cm _2 .s ) that follows the diffusion equation: 



1 drii 1 1 dT 
tij dz Hi T dz 



-mK 



Id* 

y» dz 



(5) 



where K is the eddy diffusion coefficient (cm 2 .s _1 ), D, is the 
molecular diffusion coefficient (cm 2 .s~') and //, the scale height 
of the species i. 

At both upper and lower boundaries, we impose a null flux for 
each species. 



2.4. Photochemistry 

We add to the thermochemical scheme a set of 34 photodisso- 
ciations, which are presented in Appendix [D] As we can see 
in Fig. S for HD 209458b and HD 189733b, UV flux pene- 
trates down to a pressure of about 1 bar, where the tempera- 
ture is higher than 1500 K. At these temperature and pressure, 
endothermic reactions do matter, which implies that photochem- 
istry and thermochemistry are coupled in such highly irradiated 
atmospheres. We used absorption cross section at the highest 
available temperature (i.e. 370 K at maximum, which is low 
compared to the temperatures in the atmosphere of hot Jupiters 
(see Fig.|4|i). 

To calculate the photodissociation rates in all the layers of 
the atmosphere, we compute the stellar UV flux as a func- 
tion of pressure and wavelength, taking into account molecu- 
lar absorption by 22 species (Appendix |D]i and Rayleigh scat- 
tering. Actinic fluxes are calculated with a resolution of 1 nm 
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Fig. 3: Penetration of UV flux in the atmosphere of HD 209458b 
and HD 189733b at the steady state in function of wavelength. 
Plots represent the level where the optical depth r = 1 . The name 
of the compounds responsible for the main absorption at differ- 
ent wavelengths is indicated. 



(which is also the resolution we adopted for the absorption 
cross-sections), assuming a plane parallel geometry a nd an in- 
cidence angle 9 of 48° (as done in iMoses et al. (2011), because 
< cos 9 >= 2/3 (9 - 48°) is the projected-areaweighted average 
of the cosine of the stellar zenith angle over the planetary disk 
at secondary-eclipse conditions). Multiple Rayleigh scattering is 
coupled w ith absorption throu gh a simple two-stream iterative 
algorithm (llsaksen et al.ll 1977b . 
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2.5. Application to hot Jupiters: HD 209458b and 
HD 189733b 

HD 209458b and HD 189733b are transiting planets around 
nearby bright stars. Their atmospheres have been studied 
by their transmission spectrum obtained during the pri- 
mary transit and their day-side emission spectrum mea- 
sured at the secondary eclipse. These observations can be 
used to constrain the therm al profile (Swa in et a D |2009a; 
Madhusu dhan & Seagerl [2009b and t o detect the spectral sig- 
nature of at mospheric compounds jC harbonneau et al. 2002; 
Tinetti et al.1 120071: ISwain et alj l2008ri iGrillmair et al. 2008; 



Table 1: Properties of the systems HD 209458 and HD 189733 



HD 209458 



HD 189733 



Landand-Shula et al. 2009; S wain et al . 2009a b; Bea ulieu et all 
20101) . 

In this preliminary study, we do not compare the results of 
our model with these observations for two reasons. First, be- 
cause there is not yet a consensus on the actual constraints that 
can be drawn from these measurements. Secondly, such a com- 
parison would imply to address the effects of circulation on the 
composition and to explore all the range of possible elemental 
abundances for these objects. Several recent works claim that ob- 
servati ons of some hot Jupiters im ply enhanced elemental C/O 
ratios (Madhusudhan et all 1201 lallbh . With C/O ratios close or 
above unity, species with more than two carbon atoms will be 
important and our C0-C2 network does not allow us to study 
them. For these reasons, we are implementing a C0-C6 mech- 
anism and a coupling with atmospheric circulation, which will 
be described in further studies. At this stage, our main goal is to 
compare the results of ou r mode l with already published works, 
in particular iMoses et alj (I20T1 . hereafter Mil. We could also 
have compared our results with those of IZahnle et alJ (|2p09a), 
who explored a broader range of conditions, included sulfur- 
bearing species and various elemental compositions. We decided 
to restrict our comparison with Mil because their model, like 
ours, only includes species made of C, H, O and N (and He), and 
also because Mil already did a co mparison between their results 
and those of lZahnle et alJ d2009al) showing only little discrepan- 
cies when the same conditions are considered. We used the same 
conditions (P-T profiles, eddy diffusion, elemental abundances) 
as in Mil, so that differences should come only from kinetics 
(and photochemistry in the upper atmosphere), which represents 
the novelty of our approach. 

2.5.1. Physical properties and composition 

The physical prop erties of HD 209458b have been refined by 
iRowe etalJ (l2008) and are presented in Table[T] with some prop- 
erties of th e host star. Properties of H P 189733b and HD 189733 
come from lSouthworfhl d2008l 120 1 Oh . 

In order to compare the outcomes of the two models 
(Sect. 13.21 ). we use the temperature and eddy diffusion profiles 
published in Ml 1 (Figs.[4]i. Als o following M 1 1 , we as sume pro- 
tosolar elemental abundances (Lodders & Palme 2009) for both 
planets, with 20% of depletion for oxygen (sequestered along 
with silicates and metals). We start our time-dependent model- 
ing with the thermodynamic equilibrium abundances calculated 
with TECA (an equilibrium model described in Appendix iBl) at 
each level of the atmosphere. 



2.5.2. U V spectral irradiance 

As HD 209458 is a GO star (Table [B, we use the UV spectral 
irradiance of the Sun for this star. For the star HD 189733, a 
Kl - K2 star (Table Q]), the UV spectrum has been provided 



Distance Sun-Star (pc) 


47 


19.3 


Distance Planet-Star (AU) 


0.047 


0.03142 


M p {M Jup ) 


0.69 ± 0.01 


1.150 ±0.028 




1.339 ±0.002 


1.151 ±0.036 


PpiPjupiter) 


0.26 ± 0.04 


0.755 ± 0.066 


go(m.s- 2 ) 


9.54 ± 0.69 


21.5 ± 1.2 


P orbil (days) 


3.5247489(2) 


2.21857578(80) 


MM) 


1.083 ±0.005 


0.840 ± 0.030 


R*(Re) 


1.118 ±0.002 


0.752 ± 0.023 


Spectral type 


GOV 


Kl V - K2 V 



Note.- The ltr uncertainty in P or bit is given in parentheses in units 
of the last digits 



to us by Ignasi Ribas (private communication). It is based on 
FUSE and HST observations of the star e Eridani, a proxy of 
HD 189733 (similar type, age and metallicity), in the 90-330 
nm range. Between 0.5 and 90 nm, i t is based on data from th e 
X-exoplanets Archive at the CAB dSanz-Forcada et all 1201 ll) . 
Above 330 nm, we use a synthetic s pectrum calculated wit h 
the stellar atmosphere code Phoenix (iHauschildt et al.1 11999). 
This model for the UV spectrum of HD 189733 slightly differs 
from the one chosen in Ml 1. We also tested our model with the 
spectrum used by Mil and found negligible differences at the 
pressure levels we model. 



3. Results 

3. 1 . Nominal model 

First of all, we checked that our kinetic model reproduces the 
thermodynamic equilibrium, in the absence of vertical mixing 
and photodissociation. We obtained differences lower than a few 
percent, as found with the 0D model (see Sect. I2.2I ). For both 
planets, the homopause is always found above the 1 x 10~ 5 mbar 
level, which is beyond the range of pressure that we model. As 
a consequence, and although it is included, molecular diffusion 
does not affect our results. 

Figure [5] shows the steady-state composition of the atmo- 
sphere of HD 209458b and HD 189733b, with vertical transport 
and photodissociations, while in Fig. [6] photodissociations have 
been removed. Comparing Figs. [5] and [6] shows us the influence 
of photolyses. Although HD 209458b receives a higher UV flux 
than HD 189733b, we can see that UV photons have little effect 
on the composition of HD 209458b, while they have a signifi- 
cant influence on the chemistry of HD 189733b. This is because 
the temperature is higher in HD 209458b so that the chemical 
timescales are significantly shorter than the lifetime of species 
against photolyses. So in HD 209458b, kinetics dominate over 
photodissociations, even at high altitude. In HD 189733b, how- 
ever, photodissociations affect the composition down to about 
the 10 mbar level. This is particularly noticeable for H and OH 
abundances. The production of H is dominated by the photolysis 
of H 2 for pressures lower than 1 yubar. Below this level, and for 
pressures higher than 0. 1 mbar, H is produced by the photodis- 
sociation of H 2 0, with a minor contribution of the photodisso- 
ciations of NH 3 and HCN. The abundance of OH follows the 
profile of H, and increases for pressures lower than 10 mbar. 
There is a photochemical enhancement of HCN above the 10 
mbar pressure level, as discussed in Mil. CH 4 is destroyed by 
photolyses for pressures lower than 0.01 mbar. NH 3 is photodis- 
sociated down to levels as deep as 1 bar, but vertical transport 
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Fig. 4 : Pressure-temperature profiles (left) and eddy diffusion profiles (right) of HD 189733b and HD 209458b (from lMoses et alJ 
12011b . 



compensates this destruction for pressures higher than 0. 1 mbar. 
Above that level, the amount of NH 3 decreases with altitude due 
to photolyses. Photoch emistry has a negligible effect on CO z , as 
noted bv lZahnle et all d2009al) . 

For HD 209458b, we can see in Fig. |6]that mixing quenches 
NH 3 and HCN at 1 bar and CH 4 at 400 mbar. These species are 
transported up to ~1 mbar pressure level, but as the tempera- 
ture increases with altitude at this level, they tend to come back 
to their thermochemical equilibrium values, so their abundances 
decrease again. For the other molecules, like H 2 , H and C0 2 , 
at the thermodynamic equilibrium, there is a steep variation of 
composition (smoothed by vertical mixing) corresponding to the 
strong temperature gradient of the upper atmosphere. 
Vertical quenching has an effect on a larger part of the atmo- 
sphere of HD 189733b. NH 3 and HCN are quenched at 5 bar, 
CH 4 at 1 bar, H at 40 mbar and C0 2 at 20 mbar. Quenching 
contaminates the composition up to very low-pressure levels 
(lO^mbar). 

3.2. Comparison with Moses etal. (2011) 

3.2.1. Equilibrium 

Overall, the composition we calculate at thermodynamic equi- 
librium (which is our initial condition) is very close to what is 
obtained by Ml 1 except for HCN for which there is a difference 
that can reach ~ 30% at 100 bars and 1545 K. To check our cal- 
culations we also did a comparison with the code STANJAN0 
and we found negligible discrepancies for the species we com- 
pared, including HCN at this pressure and temperature. The dif- 
ference with Mil probably comes from the coefficients used for 
the NASA polynomials (see Appendix lAl. Although the differ- 
ence remains small for equilibrium calculations, we should keep 
in mind that it may significantly affect the kinetics of HCN and 
related species through the calculation of the rate for backward 
reactions and vertical quenching. 

3.2.2. Steady-state 

In Fig. 13 we compare our results at steady state with those 
of Mil. Differences between Mil and our nominal model are 



5 http://navier.engr.colostate.edu/tools/equil.html 



also shown species by species in Figs. [8] and [9] Discrepancies 
between the two models can be due to the different chemical 
schemes and, at levels where the results are sensitive to photo- 
chemistry, to possible differences in the UV fluxes, cross sec- 
tions, and photodissociation quantum yields. An influence of the 
numerical implementation (like the discretization of the atmo- 
sphere, the solver for the continuity equations, or the treatment 
of the UV transfer) is also possible. 

In the lower atmosphere of HD 189733b and for most of 
the atmosphere of HD 209458b, photolyses have a negligible 
influence and departures should be caused by the kinetics. For 
these regions we find very similar results for species that re- 
main at their equilibrium abundance (H, OH, CO, C0 2 , H 2 
for instance), which only confirms, as stated before, that our 
thermodynamic equilibrium codes are in good agreement. For 
species quenched by mixing, however, significant deviations ap- 
pear, in particular for NH 3 , HCN and CH 4 . Their quenching oc- 
curs at different pressure levels and, thus, for different abun- 
dances that will then contaminate a large fraction of the at- 
mosphere above. The discrepancies are much more significant 
in the case of HD 189733b, due to a higher sensitivity to ki- 
netics. Although the kinetic network is certainly the main rea- 
son for these departures, it is also true that quenching can be 
quite sensitive to the resolution of the pressure (or altitude) 
grid, in particular when there is a steep gradient of tempera- 
ture which is the case in the convective zone (P >100 bar). 
For this reason, we impose the thickness of individual lay- 
ers to be smaller than 1/8' A of the local scale height, which 
results in ~300 layers for the pressure range that we model. 
Although we do not know what resolution is used in Mil, it 
seems more likely that the deviation comes from differences 
in the kinetic network itself. As explained in the Introduction, 
we use a chemical scheme validated for the species represented 
and for most of the range of temperature and pressure of the 
modeled atmospheres. Mil, on the other hand, use a chemical 
scheme derived from Jupiter and Saturn m o dels {G ladston e et al] 
1 19961 lMoses1ll996t iMoses et al.lll995b1laL l2000allbh completed 
by hi gh-temperature kinetics from combustion-chemistry litera- 
ture dBaulch et al II 19921 [l994ll2005t lAtkinson et al.lll997il 2006; 
Smit h et al.lll999UTsan|dll987l 11991b iDean & Bozzellil feoOO). 
which has not, to our knowledge, been validated against experi- 
ments. We also note departures in the upper atmosphere, where 
photolyses are important. In particular, H and OH exhibit similar 
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Fig. 6: Steady-state composition of HD 209458b (left) and HD 189733b (right) calculated with our nominal model without pho- 
todissociation (color lines), compared to the thermodynamic equilibrium (thin black lines). 




profiles than those of Ml 1 in HD 189733b but shifted by about 
one order of magnitude in abundance for pressures lower than 
50 mbar. CH 4 is also affected. We checked that these differences 
are not due to the use of different stellar fluxes by switching be- 
tween the flux we use and the one used in Ml 1 (for HD 209458b, 
we both use the solar UV flux). At these altitudes, we note a sig- 
nificant sensitivity of the mixing ratio of these species to the 
Rayleigh scattering, so the treatment of the scattering could ex- 



plain at least part of this disagreement. Again, and although we 
do not know the details of the photochemical data and radiative 
transfer used in Ml 1, we assume that kinetics explain the differ- 
ences. 
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3.3. Other networks for nitrogen species 

The main differences between Mil and our results is related 
to the quenching of NH 3 and HCN. As mentioned in Mil, the 
chemistry of nitrogen compounds has been less studied than car- 
bon species and chemical networks have been subjected to less 
validation. However, NOx, HCN, CN and NH 3 are important 
species in applied combustion (gas fuel, for instance, can con- 
tain high concentrations of ammonia), and should be well re- 
produced within the temperature and pressure range of the vali- 
dation. Quenching is found to occur within 1 to 10 mbars, cor- 
responding to the range of validation in terms of pressure. An 
originality of our network compared to other schemes used in 
combustion is that it is not optimized to increase the agreement 
between modeling and experiments. In other words, the rate co- 
efficients of the individual processes have not been altered com- 
pared to their original measurement or estimate. The application 
of the network is therefore not strictly restricted to the validation 
domain. 

Other sub-mechanisms are available to model the kinetics of 
nitrogen-bearing species. They have been constructed based 
on different approaches (optimization, specific domain of ap- 
plication, reduced number of reactions). In order to test our 
model against other nitrogen schemes, we replaced our nitro- 
gen reaction base by nitrogen sub-mechanisms taken from other 
C/H/O/N mechanisms: 

- GRIMECH mechanism based on GRI-Mech 3.0 
( Smit h et al.| [l999) with several reactions involvi ng NOx com- 
pound s added with respect to the mec hanisms of [G laude et all 
(2005) as recommended and done bv lAnderlohr et alJ (12009). 
It includes 162 reversible reactions involving 26 nitrogen com- 
pounds. The GRTMech 3.0 is a mechanism designed to model 
natural gas combustion, including NO formation and reburn 
chemistry. As already mentioned in Sect. 12.1.21 it has been 
optimized as a global mechanism, i.e., some rate coefficients 
have been modified (relatively to the literature) in order to 
fit the results of a pool of experiments with conditions and 
compositions specific to combustion. The individual processes 
have not been studied separately in all the pressure and tem- 
perature range. Applying this mechanism beyond its domain of 
optimization/validation is a risky extrapolation. Mixing ratios of 
oxidants, for instance, are very low in hot Jupiter atmospheres 
compared with the experiments used to optimize/validate 
GRI-Mech 3.0. 

While doing the present study, we noticed that two reactions 
from the GRIMECH mechanism had wrong rates and corrected 
them. These erroneous rate constants, that can be traced bac k 
to the NIST Chemical Kinetics Database dManion et al.l l2008). 
were identified by systematically comparing reaction rate 
constants with collision limit values and energy barriers with 
the enthalpy budget of the reaction. 
The first reaction is 



N 2 + H 



N( 4 S) 



N( 4 S) 



H 



(6) 



for which the rate given by NIST is k f (T)=l.26 x 10" y 
y-0.20 e -27,254/r cm 3 mo i ecu i e -i s -i (with T in Kelvin), although 
this expression corresponds in fact to the reaction [7] the ther- 
mal dissociation of N H through collisions with atomic nitrogen 
(ICaridade et al.ll2005l) : 



The above expression overestimates by many orders of magni- 
tude the rate constant of the reaction, whose activation energy 
must be around 100,000 K, as implied by the bond energy of 
molecular nitrogen and by measurements of the thermal dissoci- 
ation of N2 through collisions with various bodies. 
We finally adopted for the dissociation of N 2 a more general 
form: 



N 2 +M 



N( 4 S) + N( 4 S) + M 



(8) 



with the reaction rate constant kf{T)-\.66\ x 10 
r -3.30 e -n,3io/r cm 3 molecule 1 s 1 dThielen & Rofhl fl986h . 
This rate had a strong influence on our results. Note that using 
the wrong rate has a large effect on the atmospheric profiles of 
NH 3 and HCN. 
The second reaction is 



HONO + NO 



NO, + HNO 



(9) 



for which the reaction rate constant given by NIST (kj(T)-l .34 
xl0 -20 T 2.64 e -2034/r cm 3 mo i eC ule~ V) was in fact the rate o f 
the reverse of reaction (|9), as calculated by Meb el et aD (Q998). 
In our model, we use in fact the reaction : 



N0 2 + HNO 



HONO + NO 



(10) 



-12 



NH + N 



H + N, 



(7) 



with the reactio n rate constant k/(T)= 1.00 xlO 
molecule" V dTsang & Herronlll991l) . 
When this paper is published, these rates should be corrected in 
the NIST Chemical Kinetics Database, but one should check that 
the wrong rates are not used in modeling. 

- GDF-Kin , a mechanism optimized for natural ga s combus- 
tion modeling (Tur biez et al.lll998tlDe"Ferrieres et al]|2008l) that 
includes less individual processes: 180 reversibles reactions in- 
volving 22 nitrogen species. Several experimental data on natu- 
ral gas combustion have been acquired in partnership with Gaz 
de France to develop this mechanism. NOx che mistry has been 
included in GDF-Kin 3.0 (lEl Bakali et al.ll2006l). We use the up- 
date version GDF-Kin 5.0 (Lamoureu x et al.l 20101) . in which 5 
reactions involving NCN have been refined in order to better re- 
produce the kinetics of this species. It is validated for tempera- 
tures between 400 and 2200 K and pressures between 0.04 and 
lObars. 

- DEAN, taken from lDean & Bozzellilj gOOO). This book that 
presents a catalog of reactions is used by Moses et al.1 (1201 ll) . at 
least for some reactions. The mechanism derived from this work 
includes 370 reversible reactions involving 49 nitrogen species 
and one C/H/O species that is not included in our C0-C2 scheme: 
HCOH. The purpose of the work of Dean and Bozzelli was to 
list gas phase reactions involving nitrogen-bearing species that 
could be important for high temperature combustion modeling 
and to provide the associated rate coefficients based on an anal- 
ysis of elementary reaction data, when available, or on estima- 
tions from thermochemical kinetics principles otherwise. This 
mechanism was developed on the basis of analysis of individual 
reactions rather than by attempting to reproduce any specific set 
of experiments. It is clearly written in this book that: "Although 
we show in the chapter that this mechanism provides a reason- 
able description of some aspects of high-temperature nitrogen 
chemistry, we have not attempted a comprehensive comparison". 
Therefore, this kinetic network should be viewed as a database 
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Fig. 8: Abundances of C H,, HCN, NH ? and CH, in HD 209458b (left) and HD 189733b (right) with the four different models, 
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Fig. 9: Abunda nces of CH-, H, OH and H 2 in HD 209458b (left) and HD 1 89733b (right) with the four different models, compared 
to the results of iMoses et 
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of reaction rate constants rather than a validated mechanism, in 
the absence of validation. 

The impact of the different nitrogen sub-mechanisms on the 
abundance profiles of various species is illustrated in Figs.[8]and 
[9] Thermodynamic equilibrium is the same for all schemes (for 
the species in common). 

First, we restrict our analysis to pressures higher than ~ 
1 mbar in order not to mix effects due to kinetic rates with 
possible differences in the photochemical data or modeling. For 
HD 209458b, the main species that are significantly affected at 
these pressure levels by the change of nitrogen scheme are HCN 
and NH 3 . This is not surprising as these are the most abundant 
nitrogen species departing from equilibrium due to quenching, 
and the pressure level at which the quenching occurs depends on 
the kinetic network adopted. These two species (but also others) 
show even larger differences in the case of HD 189733b, as the 
lower temperatures of this atmosphere enhance the differences 
due to the kinetics. Departures between schemes are expected to 
become even more important for cooler atmospheres. None of 
the tested schemes shows an general improvement of the agree- 
ment with Mil. Similar NH 3 quenching is found by Mil and 
with DEAN for both p lanets, which makes sense as Mil use 
iDean & BozzeliU ((2000) as a source of reactions and associated 
rates for N-bearing compounds. This similarity is found also for 
HCN but only at pressures higher than 1 bar. At higher altitudes, 
the HCN profiles from Mil become closer to the result with 
GRIMECH. For both planets and both HCN and NH 3 , profiles 
obtained with GDF-Kin are bracketed by those from the nominal 
model and DEAN, while GRIMECH gives significantly higher 
abundances than all other models in the quenching region. With 
GRIMECH, we also notice that NH 3 becomes the main nitrogen- 
bearing species from the bottom of the atmosphere up to 0.03 
mbar because of vertical mixing, whereas thermodynamics pre- 
dicts that N 2 should be the main nitrogen-bearing species. 

Understanding the roots of theses discrepancies would 
require an in-depth study of the sensitivities of these schemes 
to reaction cycles, as a function of temperature and pressure, 
which is far beyond the scope of this study and would require 
tools that may have to be developed specifically for such large 
networks. To illustrate this difficulty let us identify the reaction 
that dominates the production rate of NH 3 for HD 189733b, at 
100 mbar. We find it to be the same for all mechanisms: 



NH 2 + H 2 



NH 3 + H, 



(ID 



whose rate constant is similar in all schemes, and is calculated 
reversing the reaction: 



NH 3 +H 



NH 2 + H 2 , 



(12) 



which dominates the destruction of NH 3 , also in all the schemes. 
The following rates are found in the different schemes: 

- Nominal, Dean, GRIMECH : 9.00 x iQ- i9 r 2.4 e -4990/r 
cm 3 . molecule -1 . s -1 , derived from lKo et al.l (1 1990b . 

- GDF-Kin: 1 056 x lQ- 18 r 2 - 3 V 5114/:r cnrlmolecule^.s" 1 from 
Mic hael et al.1(ll986h . 

We could think that these slight differences are responsible 
for the different results. However, NH 3 does not display the 
same abundance when using Dean and GRIMECH mechanisms, 
despite the fact that they both share the same rate constant. 
Moreover, nulling the rate constant of this reaction in the nom- 
inal scheme, does not affect the quenching level of NH 3 , nor 
its abundance for pressures higher than 10 mbar. We can there- 
fore eliminate this hypothesis. Key reactions are in fact usually 



those that limit the rate of a cycle and that do not dominate the 
production or destruction of a given species. Finding those lim- 
iting processes in complex networks is a field of research in it- 
self. Identifying key pathwa ys and their l i miting reactions re- 
quire dedicated algorithms dLehmannl2004t Grenfell et al . 2006; 
iDobriievic et alj|2010t fStock et all 1201 U 1201 2|) whose adapta- 
tion to the large networks we consider will require further work. 

For hydrocarbons (see for instance CH 4 , C 2 H 2 and CH 3 ) all 
the models we tested cluster to the same profiles for pressures 
below 1 mbar. This shows that N-bearing species have little in- 
fluence on hydrocarbon chemistry at these altitudes. (This would 
no longer be true at higher temperature and for higher C/O ra- 
tios as HCN would become a major reservoir of both N and C). 
Mil systematically finds higher mixing ratios for hydrocarbons 
(but within one order of magnitude) above the quenching level 
of CH 4 (1-10 bar), likely due to kinetic differences in the C0-C2 
mechanism. 

At lower pressure, Figs. [8]and|9] show large differences that 
are no longer due to quenching. At pressures lower than 1 mbar, 
the abundances of hydrocarbons depend on the nitrogen network 
used. It is particularly striking for C 2 H 2 in HD 189733b, where 
DEAN and GRIMECH, on the one hand, and the nominal model, 
GDF-Kin and Ml 1, on the other hand, seem to cluster in two dif- 
ferent regimes, exhibiting 2 to 3 order of magnitude differences 
at 0.1-0.001 mbar. Departures between network results can be 
due to differences in the kinetic network (different reactions, dif- 
ferent rates, different minor species included) but also in photo- 
chemistry. Indeed, some UV-absorbing species are not included 
in all the models, such as N 2 H 4 , HN0 3 , C 2 N 2 and N 2 4 , which 
have absorption domains that overlap that of C 2 H 2 . 



3.3.1. Corresponding emission and transmission spectra 

In order to calculate the planetary transmission and emis- 
sion+reflection spectra of HD 189733b (Fig. ITOb and 
HD 209458b, we u se a line-by-line radiative transfer m odel 
from 0.3 to 25 fim (llro et all 120051 llro & Demingll2010h . The 
opacity sources included in the model are the main molecular 
constituents: H 2 0, CO, CH 4 , NH 3 , Na, K and TiO; Collision 
Induced Absorption by H 2 and He; Rayleigh diffusion and H 
bound-free and H 2 free-free. For absorbing species not included 
in our kinetic model (Na, K and TiO), chemical equilibrium is 
assumed. The current model does not account for clouds. For 
the reflected component, we use synthetic stellar spectra gener- 
ated from ATLAS0 The main difference from the static model 
described in llro et al.l (|2005) is the additio n of NH, for which 
we used the HITRAN 2008 database dRothman et all 12009). 
Planetary parameters are taken from Table [TJ 

We applied this model for the compositions obtained with 
the two nitrogen mechanisms, which give the most opposite re- 
sults (Nominal and GRIMECH), as well as for chemical equi- 
librium. The GRIMECH scheme gives the highest abundance 
for ammonia: ten and one hundred times more NH 3 than the 
nominal model for HD 209458b and HD 189733b, respectively. 
As a consequence, features of this molecule become noticeable 
on both the emission and transmission spectra at 1.9, 2.3, 3.0, 
4.0, 6.1 and 10.5 yum. The most prominent feature is found for 
HD 189733b at 10.5 fim. NH 3 features are also visible on the 
spectra of HD 209458b, but so slightly that it would not be ob- 
servable. 

At the moment, our radiative transfer model does not include 
the contribution of HCN to the opacities. Based on the HCN 



http://kurucz.harvard.edu/stars.html 
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abundances and associated spectra found by Ml 1 , we can expect 
the spectra to be also sensitive to the HCN abundance. Indeed, 
at the altitudes probed by the observations, there is nearly two 
orders of magnitude less HCN with our nominal model than with 
GRIMECH, and GRIMECH gives HCN abundances similar to 
that of Mil. Therefore, the signature of HCN found by Mil 
at 14 fim should also become noticeable with the GRIMECH 
version of our scheme. 

Some observational data are superimposed to the spec- 
tra dCharbonneau et al.ll200& iGrillmair et al]l2008t ISwain et alJ 
2009 b|) (for emission spectrum) and dKnutson et alj 



2007 



2009; Swain et al. 2008b; Beaulieu et al. 2008; Sing et al. 2009 
Desert et aLll2009t lAgol et alJl2010h (for transmission spectra), 
but we do not discuss the agreement between observations and 
synthetic spectra, as we did not attempt to fit the observable us- 
ing different thermal profiles. 

4. Discussion 

To study HD 209458b and HD 189733b, we use the 
pressure-temperature and eddy diffusion profiles from Mil. 
These profiles are der ived from general circulation models of 
Showman et a D (11009). This choice is motivated by the com- 
parison with Mil and also because the actual physical struc- 
ture of these atmospheres is not yet well constra ined by obser- 
vations. Note however, that iHuitson et alJ (1201 2l) recently pub- 
lished thermal profiles for both planets, as inferred from tran- 
sit spectroscopy. It seems that HD 189733b could be warmer 
than HD 209458b between 10~ 5 and 10 3 bar. Considering the 
large uncertainties affecting their physical conditions, our mod- 
els should not be considered as predictions of the composition 
of hot Jupiters but more as a step in developing chemical models 
of these objects and model intercomparison. 
In addition, circulation is not included in our model nor, to 
our knowledge, in other current photochemical models of hot 
Jupiters although it has a significant influence on the chemical 
composition of the atmospheres due, for instance, to the strong 
longitude dependency of the temperature. We study the influence 
of the horizontal transport on the composition of the atmosphere 
of HD 209458b in a forthcoming paper (Agundez et al., in prepa- 
ration). 

Because modeled abundances depend significantly on the re- 
action network (in particular for NH 3 , HCN and some hydrocar- 
bons), we recommend using a network that has been validated 
(but not optimized) against experiments for conditions as close 
as possible to those of application. We do not claim that the net- 
work we release is a definitive one, it will necessarily evolve 
as new experimental results and kinetic/thermodynamic data be- 
come available. Detailed nitrogen chemistry, for instance, has 
been implemented in combustion networks only recently and 
will be subjected to further evolution. Missing elements should 
also be added to the network, sulfur being the most obvious one. 
Since the scheme has not been optimized, adding new species 
and reactions to the network we release is possible. Moreover, al- 
though the range of P, T and the elements considered in combus- 
tion model do fit well with the study of hot Jupiter atmospheres, 
the ratio between hydrogen and the other elements does not. For 
this reason, among others, it is important to avoid using opti- 
mized networks that would prevent modeling of such hydrogen- 
rich mixtures. Note also that our current network, that cannot be 
used to study the abundance of species with more than 2 carbon 
atoms, is likely insufficient to study atmospheres with C/O ratios 
close to or above unity. For this reason we are currently working 
on an extended network that can model species up to 6 carbon 



atoms. 

Hot Jupiter atmospheres represent an extreme case of planetary 
atmospheres in terms of both high temperatures and low metal- 
licity. With the progress of observations, cooler and heavier at- 
mospheres are or will be soon (with JWST, EChO, Finesse, E- 
ELT) accessible to characterization. Cooler atmospheres will de- 
part more from equilibrium and will thus be more sensitive to the 
details of the kinetics. Molecules that remain minor constituents 
of hot Jupiter atmospheres will become more abundant in cooler 
atmospheres and have an increased influence on their spectral 
appearance and thermal structure, making the use of validated 
schemes even more relevant. More metallic atmospheres should 
be found as we explore the exoplanet realm towards smaller ob- 
jects with higher core/envelope masses, and eventually terrestrial 
objects. While some uncertainties still exist when applying our 
network to hydrogen-rich atmospheres, atmospheres with de- 
creasing hydrogen to heavy elements ratio become closer to the 
conditions of validation (by their equivalence ratios, to use com- 
bustion terms). Even hot Jupiters can be significantly enriched 
(by factor of 10 or more) in heavy elements compared to their 
parent stars, which c ould for instance ex plain the high observed 
abundances of CO z dZahnle et al.ll2009ah . 

5. Conclusion 

We have constructed a chemical scheme to study the atmo- 
spheric composition of hot Jupiters. Compared to existing mod- 
els we chose to use one that 

- is derived from mechanisms that are intensively used for indus- 
trial applications (in particular car engine simulations), 

- has been subjected to in-depth validation protocols in a broad 
range of temperatures, pressures and compositions, 

- is based on individual rate coefficient that have not been altered 
in order to optimize the agreement between the collective be- 
havior of the network and experiments (contrary to most mech- 
anisms in combustion). This allows users to apply the network 
slightly beyond its validation domain and to add species and re- 
actions, 

- uses experimental measurements for some of the endothermic 
reactions when robust data are available but still reproduces ther- 
modynamic equilibrium with an excellent accuracy. 

We developed a ID photochemical model based on this 
kinetic scheme, and which includes vertical transport (mixing 
and molecular diffusion) and photodissociations. We applied 
this model to the hot Jupiters HD 20945 8b and HD 189733 b 
and compared our results with those of Mos es etaD d201 ll) . 
Qualitatively, we find similar conclusions: photodissociations do 
not have a significant impact on the atmospheric composition of 
HD 209458b, with the high temperatures we assume. It remains 
at the thermodynamic equilibrium for pressures higher than 1 
bar. For lower pressures vertical transport affects the abundances 
of HCN, NH 3 and CH 4 and some of the minor species associ- 
ated. For HD 189733b, we assume significantly lower tempera- 
tures and find the atmospheric composition to be more sensitive 
to photolyses and vertical transport, all species being affected, 
except the main reservoirs, H 2 , H 2 0, CO and N 2 . Quantitatively, 
however, we find significant differences (up to several orders of 
magnitude in the case of HD 189733b) in the abundances that are 
likely to be due to the different chemical schemes used. These 
differences are smaller for HD 209548b because kinetics have 
less influence. The quenching of HCN and NH 3 , as well as CH 4 
to a lower extent, is particularly affected, as well as most species 
sensitive to photochemistry in the upper atmosphere. Despite be- 
ing large in terms of abundances, these differences do not pro- 
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Fig. 10: Synthetic day-side (left) and transmission (right) spectra of HD 189733b with the nominal mechanism (green curve) com- 
pared to the one corresponding to the GRIMECH mechanism (red curve) and to the thermochemical equilibrium (blue curve). 
The dark curve is obtained when NH 3 is removed from the model. The day-side fluxes are given as brightness temperatures (T/,). 
Because of the reflection component, note that the link between T/, and the atmospheric thermal profile is altered below 2 fim. The 
transmission spectrum is given as the apparent planetary radius. The data points obtained from various observations are also shown. 



duce strong effects on the spectra, with the exception of NH 3 
at 10.5 pm. Confronting different schemes with observations 
will thus have to wait more accurate spectroscopic observations 
(JWST, E-ELT, EChO). Until that, experimental validation ap- 
pears mandatory. 

In order to illustrate the sensitivity to the kinetic scheme, we im- 
plemented different available nitrogen schemes, that are either 
optimized (GRIMECH, GDFKin) or non validated (Dean). We 
studied the extent of the possible results, and found large differ- 
ences whenever disequilibrium chemistry is at work. Changing 
the nitrogen scheme strongly affects the quenched species 
(HCN, NH 3 ) and most species (including hydrocarbons) in the 
upper atmosphere of HD 189733b. For HD 209458b, deviations 
are again less noticeable as the atmosphere departs less from 
equilibrium. We therefore emphasize on the need to use vali- 
dated and non optimized chemical schemes. This is already true 
for hot Jupiters but this is even more crucial in the case of cooler 
atmospheres (GJ1614b, GJ3470b, for instance), which depart 
more from thermodynamic equilibrium and are more sensitive 
to kinetics. 

Our nomi nal scheme can be downloaded from the KIDA 
database^ ( Wakela m et al.ll20T2l) . The scheme is designed to re- 
produce the kinetic evolution of species with less than two car- 
bon atoms. In order to study atmospheres with C/O ratio higher 
than solar (close to or above 1), we are currently developing a 
Co-C6 scheme, which will be able to describe kinetics of species 
up to 6 carbon atoms. One of the next improvement of our model 
should be the addition of sulfur. As kinetics of nitrogen species 
is an active field of research, we expect regular updates of the 
network (which would be notified and available on KIDA). 

Note also that conclusions of this study on the chemical com- 
position of hot Jupiters, which derive from models using an av- 
erage ID vertical profile, will probably have to be revisited with 
the effects of atmospheric circulation. 
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Appendix A: Thermochemical data 

Thermochemical properties, such as enthalpies of formation, en- 
tropies and heat capacities are of great importance to ensure the 
consistency between the rate parameters of the forward and re- 
verse elementary reactions. They are also useful to estimate the 
heat release rate. Thermochemical data for all molecules or rad- 
icals have been estimated an d stored as 14 NASA polynomial 
coefficients, according to the McBrid e" et al.l d 19931) formalism. 
The NASA polynomials have the following form: 



c° P (T) 
RT 
h°(T) 



(A.l) 
(A.2) 



T t z t 1 r 4 a b 

fll + fl7 H a-> h 04 h As 1 

RT 2 3 4 5 T 

s®(T) 

= a\ In T + a^T + ^3— + CI4— + ^5— + cli (A. 3) 

R 2 3 4 
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Where a,, i e [1,7], are the numerical NASA coeffi- 
cients for the fourth-order polynomial. Each species is charac- 
terized by fourteen numbers. The first seven numbers are for the 
high-temperature range, generally from 1000 to 5000 K, and 
the following seven numbers are the coefficients for the low- 
temperature range, generally from 300 to 1000 K. When these 
parameters are not available in the lite rature (Mc Brid e et alJ 
ll993UAllendorl2006tlGoose et alj[2010h which is the most fre- 
quent case for species present in automotive fuels, they have to 
be estimated. In this case, these data were automa tically cal- 
culated using the software THERGAS dMuller et alJI 19951) . de- 
veloped in the LRGP laboratory, based on the g roup and bond 
additivity methods proposed by Benson (1 19761) and updated 
based on the data of Domalsk i & Hearing! ([1996). The enthalpies 



of formation of alkyl radicals have been also updated accord- 
ing to the values of bon d disso ciation ener gies published by 
iTsang & Hampsonl d 19861) and bvlLuo l (|2003[) and following the 
recommendations of iBenson & Co hen ( 1997) respectively. 

An elementary reversible reaction i involving L chemical 
species can be represented in the following general form: 



v uXi 



(A.4) 



z=i 



where v' n are the forward stoichiometric coefficients, and VI are 
the reverse ones, xi is the chemical symbol of the l ,h species. 
The kinetic data associated to each reaction are expressed 

Eg 

with a modified Arrhenius law k(T) = A x T" exp «t where T 
is the temperature, E a is the activation energy of the reaction, A 
the pre-exponential factor and n a coefficient which allows the 
temperature dependence of the pre-exponential factor. If the rate 
constant associated to the forward reaction is kfj(T), then the one 
associated to the reverse reaction is k r j(T), verifying: 



K„ 



kfi(T) lk B T 



k,i(T) \ P° 



(A.5) 



where K p j is the equilibrium constant, when the activ ity of the 
reactants is expressed in pressure units dBenso 



K pi = exp 



~R~ ~ ~RT~ 



(A.6) 



where AS*' and AH9 are the variation of entropy and enthalpy 
occurring when passing from reactants to products in the reac- 
tion i, P° is the standard pressure (P° =1,01325 bar), ks is the 
Boltzmann's constant and v/ are the stoichiometric coefficients 
of the L species involved in reaction i: v; = V!- - V,.. Combined 

o A//" 

with Eqs. ( IA.2t and ( IA.3l l, -jj- and -=f can be calculated with 
the NASA coefficients: 



Appendix B: Chemical equilibrium calculation 

To compute the equilibrium abundance of the species in a def- 
inite system considered as an ideal gas, we have developed a 
Thermodynamical Equilibrium Calculator TECA. TECA is a 
software which allows the equilibrium calculation for a com- 
plex mixture. More specifically, for a given initial state of an 
ideal-gas mixture, the chemical-equilibrium program is able to 
determine the gas composition at a defined temperature and 
pressure. This calculation is based on the principle of the 
minimization of Gibbs energy (e.g. iGibbsl Il873t IWhite et all 
1958HEriksson & Rosed 1 97 UlSmith & Misserill982tlRevnoldsl 



1986): 



(B.l) 



where L is the total number of species, gi is the partial free en- 
ergy of the species / and Ni is the number of moles of the species 
I. 

The partial free energy of a compound /, behaving as an ideal 
gas, is given by: 



g} = g,(T,P) + RT\aN, 



(B.2) 



where gi(T, P) is the free energy of the species I at the tem- 
perature T and the pressure P of the system and R is the ideal 
gas constant. 

For an ideal gas, gi{T, P) is given by: 

g,(T,P) = h%T) - T$(T) + RT In (^) (B.3) 

where h^iT) and s^iT) are respectively, the standard-state en- 
thalpy and entropy of the species I at the temperature T of the 
system. 

The enthalpy and the entropy are expressed as NASA poly- 
nomials as described above. 



Appendix C: Pressure-dependent reactions 

Under some conditions, several reactions do not have the same 
rate constant depending if they occur under low or high pres- 
sure (respectively ko(T) and k^T)). In this case, between these 
two limits it appears what is called a fall-off zone. This is typ- 
ically the case in reactions requiring a collisional body to pro- 
ceed, such as thermal dissociation or recombination (three-body) 
reactions. In the present kinetic model we have different types 
of reactions with pressure dependent rate constants (Table ICTT i. 
In some cases, some species act more efficiently as collisional 
bodies than do others. Then, when available from literature, col- 
lisional efficiencies are used to specify the increased efficiency 
of the I th species in the i th reaction (see for example reaction (2) 
in TablelCHV 

For the pressure-dependent reactions, the rate constant at any 
pressure is taken to be: 



AS? ^ s°(T) 



—L = V 
R 4-i 



V,- 



and ^LV 
RT 



L U Q 
Vl 



h?(T) 



RT j-J RT 



(A.7) 



Finally, we can calculate the reverse reaction rate for the reaction 
i: 



. m k fi (T) lk B Tf-^ 



(A.8) 



k(T) = koo(T) 



P, 



l+P r 

where the reduced pressure P r is given by: 
p _ [M]ko 

r ~ l 



(CI) 



(C2) 



and [M] is the concentration of the mixture, weighted by the 
efficiency of each compound, a/, in the reaction studied: 
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[M]=£a,[X,] (C.3) 

1=1 

where [X{\ is the concentration of the species k. 

As shown in Table IC.U three methods of representation of 
the rate expression in the fall-off region are used (enhanced col- 
lisional body efficiencies of certain species are presented below 
the reaction): 

- the lLindemann et all d 19221) formulation, illustrated by reac- 
tio n (1) in Tabl e ICTl 

- the lTrod JT983) formulation, see for example reaction (2) in 
Table O 

- the SRI formulation proposed by St ewart et alJ (Il989t) . illus- 
trated by reaction (3) in Table IC.ll 



In the Lindenman form, F is unity (F =1). 
In the Troe form F is given by 



1 J7 l0 gl()(^>») 

lo Sio F = — — — (C- 4 ) 



, r logio^o+c i 

T lN-d(\og l0 (P,.)+c)\ 

with 



c = -0.4 - 0.67 x log 10 (F cent ) 
N = 0.75- 1.27 xlog w (F cent ) 
d = 0.U 



Fcent = (1 -a)exp^-^ 7 j + aexp^-| 7 j + exp|-^ r j (C.5) 

The four parameters a, T***, T* and T** must be specified 
but it is often the case that the parameter T** is not used by lack 
of data. 

The approach taken at the Stanford Research Institute (SRI) 
by IStewart et alJ (1 19891) is in many ways similar to that taken 
by Troe, but the blending function F is approximated differently. 
Here, F is given by 



F - d 



-b 



-T 



a exp h exp — 

T c 



(C6) 



where 



X - 



1 



l+(log 10 P,.) 2 



(C7) 



Appendix D: Photodissociations 
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k = AT n exp(-E/RT) 



High pressure limit 
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Table C.l: Some example of reactions with pressure-dependent rate constants present in the kinetic model. 



Pathways 




Cross sections 


Quantum yields 


Jl 


H 2 + hv 


— > H + OH 


Chan et al. (1993b1 and 


Huebner etal. (1992) 


J2 




— > H, + O('D) 


Fillion et al. (2004) and 




J3 




— > H + H + 0( 3 P) 


Mota et al. (2005) 




J4 


C0 2 + hv 


— > CO + 0( 3 P) 


Stark et al. (2007): Itvakov et al. (2008). and 


Huebner etal. (1992) 


J5 




— > CO + O('D) 


Huestis & Berkowitz (2010) 




J6 


H 2 CO + hv 


— > H, + CO 


Coooeretal. (1996) 


Huebner etal. (1992) 


J7 




— > H + HCO 


and Meller & Moortaat (2000) 




J8 


OH + hv 


— > H + 0('D) 


Huebner et al. (1992) 


Van Dishoeck & Dalearno (1984) 


J 9 


CO + ftv 


— > C + 0( 3 P) 


Olnevetal. (1997) 


Huebner etal. (1992) 


J10 


H 2 + /jv 


— > H + H 


Samson & Haddad (1994); Chan et al. (1992) 


Estimation 








Olnev et al. (1997) 




J 1 1 


CH, + hv 

4 


— > CH 3 + H 


Auetal. (1993) and 


Gans etal. (2011) 


J12 




— > 'CH 2 + H 2 


Lee etal. (2001) and 




J13 




— > 3 CH 2 + H + H 


Kameta et al. (2002) and 




J14 




— > CH + H 2 + H 


Chen & Wu (2004) 




J 15 


CH 3 + 


— > 'CH 2 + H 


Khamaganov et al. (2007) 


Parkes etal. (1973) 


J16 


C 2 H, + hv 


— > C,H + H 


CooDer et al. (1995a): Wu et al. (2001) 


Okabe (1981. 1983) 


J17 


C 2 H 3 + hv 


— > C,H 9 + H 


Fahretal. (1998) 


Fahretal. (1998) 


J18 


C 2 H 4 + fry 


> C2H2 "I" H2 


Cooper etal. (1995b) 


Holland etal. (1997) 


J19 




— > C 2 H 2 + H + H 


Orkin et al. (1997) and Wu et al. (2004) 


Chang etal. (1998) 


J20 


C 2 H 6 + hv 


— * C 2 H 4 + H 2 


Au et al. (1993): Kameta et al. (1996) 


Akimoto et al. (1965) 


J21 




— > C,H 4 + H + H 


Lee etal. (2001) 


Hampson Jr & McNesby (1965) 


J22 




— > C^H 2 + H 2 + H, 


Chen & Wu (2004) 


Lias etal. (1970) 


J23 




— > CH 4 + 'CH 2 




Mount & Moos (1978) 


J24 




— > CH 3 + CH 3 




Mount & Moos (1978) 


J25 


N, + /iv 


— > N( 2 D) + N( 4 S) 


Samson & Cairns (1964); Huffman (1969) 


Estimation 








Stark et al. (1992): Chan et al. (1993c) 




J26 


HCN + hv 


— > CN + H 


Lee (1980) 


Lee (1980) 


J27 


NH 3 + hv 


— > NH, + H 


Burton et al. (1993): Chen et al. (1999): Chene et al. (2006) 


McNesbv etal. (1962) 


J28 


NO + fry 


— > N( 4 S) + 0( 3 P) 


Iida et al. (1986): Chan et al. (1993a) 


Huebner etal. (1992) 


J29 


N,H 4 + /iv 


— > N,H 3 + H 


Vaghiiani (1993) 


Vaghiiani (1993. 1995) 


J30 


HN0 3 + ftv 


— > NO, + OH 


Sander etal. (2011) 


Estimation 


J31 


C 2 N 2 + hv 


— > CN + CN 


Benilan et al. (2012), in preparation 


Codv et al. (1977): Jackson & Halpern (1979): Eng et al. (1996) 


J32 


N,0~ 4 + hv 


— > N0 2 + N0 2 


Vandaele et al. (1998): Merienne et al. (1997) 


Sander etal. (2011) 


J33 


N^O, + /?v 


— > N0 2 + NO 


Stockwell & Calvert (1978) 


Sander etal. (2011) 


J34 


HCO + ftv 


— > H + CO 


Hochanadel et al. (1980); Loison et al. (1991) 


Estimation 



Table D.l: Photodissociations scheme used in the model 



